Bogoliubov theory of entanglement in a Bose-Einstein condensate 



Anders S0ndberg S0rensen 
Institute of Physics and Astronomy, University of Aarhus, 
DK-8000 Aarhus C, Denmark 



O 
O 

(N 

> 
O 



> 

m 
o 



I 

o 
o 



We consider a Bose-Einstein condensate which is illuminated by a short resonant light pulse that 
coherently couples two internal states of the atoms. We show that the subsequent time evolution 
prepares the atoms in an interesting entangled state called a spin squeezed state. This evolution 
is analysed in detail by developing a Bogoliubov theory which describes the entanglement of the 
atoms. Our calculation is a consistent expansion in l/\/iV, where A'' is the number of particles in 
the condensate, and our theory predict that it is possible to produce spin squeezing by at least a 
factor of l/\fN. Within the Bogoliubov approximation this result is independent of temperature. 



PACS numbers: 03.65.Ud, 03.75.Fi, 42.50.Dv 



I. INTRODUCTION 

One of the main problems in the experimental ex- 
plorations of entanglement is the engineering of a con- 
trolled interaction between individual quantum particles. 
Such an interaction among the particles is automatically 
resent in Bose-Einstein condensed dilute atomic gases 
|[ H], where the experimental results show clear sig- 
natures of the collisions of the particles. It is therefore 
natural to try to exploit this interaction to create entan- 
glement and a number of proposals have appeared which 
use the collisional interaction to entangle the atoms in 
the condensate § |, |, |, ||, |l^, 

In Ref. Q we proposed a realistic experiment where 
the atoms in a Bose-Einstein condensate are illuminated 
by a single resonant pulse. After this pulse the system 
evolves freely, and the collisional interaction among the 
particles prepares the atoms in an entangled state. The 
created entangled state even has practical applications 
in atomic clocks, where it can be used to increase the 
precision significantly In |^ we presented an ap- 

proximate simulation of the condensate dynamics, where 
the state vector was expanded on subspaces containing 
different number of atoms, and where the spatial wave- 
function in each subspace was evolved with the Gross- 
Pitaevskii equation. This calculation indicated that a 
substantial entanglement can be created, but the validity 
of the ansatz used in the calculation and the behaviour 
of the proposal at non-zero temperature is difficult to de- 
termine from the simulation. An exact calculation with 
a stochastic method ||] has later confirmed that the pro- 
posal is indeed capable of producing substantial entan- 
glement. But due to the noise in this calculation, it is 
difficult to estimate the exact amount of entanglement 
from this procedure and the method is also difficult to ap- 
ply in the regime where most experiments are currently 
operating. 

In this paper we investigate the proposal Q] in more 
detail. We first develop a Bogoliubov theory for a two 
component condensate and then we apply it to describe 
the entanglement of the atoms in the condensate. Our 
theory is a consistent expansion in the ratio between 



non-condensed and condensed particles and the validity 
of the approximation can directly be investigated for a 
given experimental configuration. To lowest order, the 
Bogoliubov theory can be used to describe the effect of 
a non-vanishing temperature, and our theory indicates 
that within the Bogoliubov approximation the created 
entanglement is independent of temperature. The re- 
sults of the Bogoliubov theory agrees very well with the 
results of the approximate method used in Ref. Q and 
thus confirms the validity of this method. 



A. Considered experimental setup 

We consider a collection of N atoms with two inter- 
nal states a and h. We assume that we first create a 
condensate consisting only of atoms in the internal state 
a. After the preparation of the condensate a fast res- 
onant pulse is applied to the atoms. The pulse is as- 
sumed to be much faster than any other time scale in the 
problem, so that the only action of the pulse is to mix 
the internal states \a) — > cos((/)/2)|a) + sin(0/2)|5) and 
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b) — sin(0/2)|a). The theory developed in 
rv| applies for all values of the angle <p but 
we only go into details with the case where the resonant 
interaction is a 7r/2 pulse, (p = 7r/2. In this paper we 
shall describe our atoms in the second quantization us- 
ing field operators and '^b, and we calculate the time 
evolution in the Heisenberg picture, where the action of 
the pulse is given by 

^a{r, t = Q+)= cos ( I ) ^'a(r , t = 0") 



(1) 



sin( ^ )*b(f,t = 0- 



*6(f , t = 0+) = cos j -^bir, t = 0-) 

+ sin(^^^*,(r,t = 0-) 

(t — 0~ and t = 0+ denotes the time just before and 
after the pulse) . After the pulse we let the system evolve 
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and perform a measurement on it at a later time. 

The resonant pulse splits the condensate into an a and 
a b component with a well defined phase between them; if 
one applies a new pulse shortly after the first, it is possi- 
ble to produce interferences between the two component 
and, e.g., transfer all the atoms back into the a state. As 
the two components evolve with time there is nothing to 
maintain the relative phase between them, and the rel- 
ative phase will spread in time because each of the two 
components do not have a well defined number of par- 
ticles ||. Due to the spreading of the relative 
phase, the magnitude of interference terms like (5']j^f,) 
will decrease with time and it will no longer be possible 
to see interferences between the two components. This 
process is referred to as phase collapse. 

Instead of focusing on the interferences between the 
condensates it is fruitful to consider other observables of 
the system. By considering fluctuations in the particle 
numbers one can show that the time evolution prepares 
the atoms in an entangled state called a spin squeezed 
stated. 



B. Spin squeezing 

Suppose that we prepare N atoms in an equal super- 
position of two internal state (|a) -I- \b))/\/2 and measure 
the number of particles in the a state. Because the parti- 
cles are independent of each other the result of the mea- 
surement will be a distributed according to the binomial 
distribution with a mean value Na = N/2 and a variance 
N/A. A spin squeezed state is a state where the atoms 
are prepared in a suitable entangled state such that the 
variance (or noise) of this measurement is reduced. 

A more convenient way to represent the atoms is to 
consider each atom as a spin- 1/2 particle with the inter- 
nal states I a) and |6) representing the spin up and spin 
down states respectively. In this language the noise in 
the counting statistics is represented by the variance of 
the Jz operator (AJ^)^, where the collective Jx, Jy, and 
Jz operators are obtained by summing the spin operators 
for the individual atoms. Expressed in terms of the spin 
operators, spin squeezing is the reduction of the noise of 
Jz (or any other spin component). Recently the first ex- 
perimental realizations of such spin squeezed states have 
been achieved 



of 



1(1, 17, m m 



In Ramsey spectroscopy, as used for instance in atomic 
clocks, a signal is recorded which is proportional to the 
length of the spin in a given direction, say the x direction, 
and the noise of the signal is proportional to the noise of 
a spin component Jg = cos(0) -f- sin(0) perpendicu- 
lar to the mean spin. The precision of atomic clocks is 
currently limited by the spin noise {AJg)^ |^, and by 
using a spin squeezed state where this noise is smaller 
than in the unsqueezed states, which are used today, it is 
possible to increase the precision of the clock. Wineland 
et al. fl^ have analysed this possibility and have shown 
that the frequency variance can be reduced by a factor 



(2) 



by preparing the atoms in a spin squeezed state with 
< 1. This points to an interesting application of spin 
squeezed states in atomic clocks. 

It has also been shown that a reduction of the squeez- 
ing parameter below unity < 1 requires the atoms to 
be in an entangled state [Q, and to quantify the entan- 
glement which may be obtained in the considered exper- 
imental setup, we shall determine the reduction of the 
squeezing parameter ^g. A different characterization of 
the squeezing is presented in Ref. pl[ |, where the depth 
of entanglement in a collection of atoms is identified. 

To see that the proposed experimental setup produces 
squeezing we shall first consider the short time evolution. 
At time t = 0+ we have £,g = 1 for all 6, and to prove 
that squeezing is produced, it is sufficient to show that 
the time derivative of £,g is negative for a suitable choice 
of 6*. 

In the limit of very low temperatures a collection of 
atoms with two internal states a and b are described by 
the second quantized Hamiltonian 



j=a,b 

+ ^ E 9n J d'r¥^ (r )^] {r)^j {r)^j (r ) (3) 



j=a,fc 

+ gab / d^ri>i{r)-^l{f)i>a{r)i'b{r), 



where Hq.j is the one particle Hamiltonian for atoms 
in state j including the kinetic energy and the external 
trapping potential Vj{r), and gjk — Ai:h^as^jk/M is the 
strength of the interaction between particles of type j 
and k, expressed in terms of the scattering length agjk 
and the atomic mass M . Here we have assumed that 
there are no spin changing collisions. A specific exper- 
imental setup where it is possible to exclude such spin 
changing collision with sodium atoms in an optical trap 
is described in 

With the Hamiltonian (^ we may find the equation of 
motion for the field operator by taking the commu- 
tator with the Hamiltonian (using h—\) 



gaa'i'l'i'a'^a 



gab^l^b^a. (4) 



For brevity most of the equations in this papers will only 
be presented for the field operator for the a component. 
The corresponding equations for the 6 component can be 
achieved by exchanging subscripts a and b. Since we work 
in the Heisenberg picture, essentially all operators and 
functions appearing in this paper will depend on time, 
and for brevity we shall in most cases omit the time 
argument on operators and functions. Unless specified 
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otherwise operators and functions, hke the field opera- 
tors and 'i'b in Eq. should always be considered 
as functions of the time t. From the derivative in (||) we 
calculate the time derivative of the squeezing parameter 
at the time t — 0^ , and by using the relation with 
= 7r/2 we can express the derivative in terms of the 
field operators before the pulse 



d 
di 



e,2=sin(20) 



gaa + 9bb — '^Qah 

27V 



d\{¥^{r,t^Q-fi>a{?,t^O-f). 



(5) 



If 9aa + dbb 7^ '^'Sab wc Can always choose the angle 9 such 
that the expression on the right hand side is negative and 
this exact equation immediately proves that squeezing 
will be produced. The time derivative is proportional to 
the two point correlation function which in the limit of 
an ideal condensate is N{N — 1) / (i^r|$|^, where $ is the 
condensate wavefunction before the pulse. 

The above argument shows that squeezing will be pro- 
duced but it does not quantify the amount of squeezing 
which is produced. A simple estimate of the squeezing 
can be obtained by making a single mode approximation 
to each of the components. In this approximation the 
part of the Hamiltonian which is responsible for spin 
squeezing can be written as ||] 



-^spin \Jz 



(6) 



with the coupling constant x given by 



^z^ gaam^ + 9bb\^T\' " 2g,fc|<i>-p|$g--p 



(7) 



where i^"}^) is the condensate wavefunction for par- 
ticles of type a (&) (the "ex" is used to indicate the exact 
wavefunction because it will later be important to distin- 
guish it from low order approximations to it). 

The squeezing arising from the Hamiltonian can 
be solved analytically [E2|. The obtainable squeezing is 
approximately « N^^, and this indicate that strong 
squeezing can be produced if the condensate contains a 
large number of atoms. The validity of the arguments 
leading to the Hamiltonian are not quite clear from 
the derivation, and it is the main purpose of this paper 
to investigate the validity of this simple description. As 
we shall see below, the Hamiltonian (|6|) is indeed a good 
approximation to the dynamics of the system but the 
coupling constant (0) has to be modified slightly. To 
obtain this result we describe the system by Bogoliubov 
theory. 



C. Bogoliubov theory 

The Hamiltonian (H) and the resulting Eq. (Q) are far 
to complicated to be solved in general and approxima- 
tions have to be applied. Here we shall determine the 



time evolution from the Hamiltonian (H) by applying the 
Bogoliubov approximation. The standard formulation of 
the Bogoliubov approximation assumes a broken sym- 
metry such that the field operators have a non- vanishing 
mean value (5'a(r,t)) = •\/iVa<I>™(r , t), where Na is the 
mean number of particles of type a. The Bogoliubov 
method proceeds by splitting the field operators accord- 
ing to 



^a(r , t) = ^aK"{r, t) + (5#„(r , t) 
^b{r,t) = y/N'b<^'irif,t) + Si'b{r,t). 



(8) 



The idea behind this splitting is that the fluctuations 5'^ a 
and S'^b are much smaller than the contribution of the 
condensates, and simple equations for the field operators 
are obtained by expanding the Hamiltonian to low 
order in the ^^'s. 

The splitting (||) gives some nice properties, such as 
simple commutation relations for the operators, but 
unfortunately it is not well suited for our present purpose. 
As we shall see below, the phase collapse and entangle- 
ment mainly arise from the evolution of the atoms in the 
condensate modes. In Eq. (||) the quantum fluctuations 
of the condensate modes and fluctuations perpendicular 
to these modes are treated equally, and this is not a good 
approximation to the system. In a direct numerical in- 
tegration we have found that the (5^s grow very quickly 
due to the dynamics in the condensate modes so that 
the assumptions of the theory break down after a short 
period of time Q . 

To circumvent the problems associated with the split- 
ting (H) we employ a different splitting which enables 
a better description of the evolution in the condensate 
modes 



(9) 



where the operator d and b obey bosonic commutation 
relations [d, d^] = 1 and [d, b] = 0. In this paper we shall 
derive the time evolution of the field operators by assum- 
ing that the fluctuations S^a and S^b are much smaller 
than the contributions from the condensate modes. An 
approach similar to ours has previously been considered 
in Refs. |2^, ^ where particle number conserving 
Bogoliubov approximations are derived. We shall gen- 
eralize these previous results to a situation where the 
atoms can be in two different internal states and where 
the number of atoms in each of the states may not be 
well defined. This generalization enables us to describe 
the collapse of the relative phase between the two con- 
densates and the entanglement of the atoms, which are 
shown to arise from the time evolution of the operators 
d and b. 

In section |ll| we present the assumptions u sed in the 
Bogoliubov description M and in section HI we de- 



rive the equations of motion for the relevant operators. In 
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section |^ we solve these equation and in ^ we apply our 
theory to the considered experimental setup described in 
subsection I A. Finally, we end the paper with a conclu- 
sion in section [v| . 



II. ASSUMPTIONS 

In this section we describe the assumptions of our the- 

Thc as- 



ory, which is presented in sections III and IV 



sumptions and methods used here are very similar to 
the two theories derived for a one component condensate 
by Castin and Dum and by Gardiner |Q. These 
two derivations are equivalent but use two different ap- 
proaches. Castin and Dum derive equations of motion 
for the relevant field operators and these equations are 
approximated by inserting a splitting as in Eq. (|) and 
neglecting high order terms in (5^'. In the work of Gar- 
diner the approximations are made to the Hamiltonian, 
where terms of high order in S'^ are neglected before the 
equations of motion are derived. Here we choose to follow 
the derivation and notation of Castin and Dum in Ref. 
[ p5[ very closely. In order to make the paper self con- 
tained we reproduce most of the calculations of Castin 
and Dum. 

We assume that the atoms are cooled to a very low 
temperature such that a condensate is formed before the 
pulse which mixes the two internal states. The conden- 
sate wavefunction is determined from the one-body den- 
sity operator which is given by 



{r'\Pas\r) = (*I(r-)^',(f' ))• 



(10) 



Mathematically, the definition of a condensate is that 
Pa.i has an eigenvalue 



(11) 



which is much larger than all other eigenvalues {N^^ ^ 
N, where N is the total number of atoms), and the spatial 
modefunction is given by ). 

The resonant pulse mixes the field operators as de- 
scribed by Eq. (^, and by using that the initial state is 
the vacuum of — 0~) we find that Pa,i{t = 0^) = 
pb,i{t = 0+) ~ Pa,i{t = 0^)/2. After the pulse the state 
|<I>°''(t = 0~)) is still an eigenvalue of the density opera- 
tor pa,i with eigenvalue 7V°'^(i = 0+) = Nl^'it = 0-)/2, 
so that we have a component of the condensate in both 
the a and the b state with wavcfunctions 

^T{r,t^ 0+) = <^T{r,t = 0+) = K^{r,t^ 0"). (12) 



and can be 



The condensate modefunctions 
used to define a splitting of the field operator as in 
Eq. (||), and to describe the squeezing we need to de- 
termine the time evolution of 5'^ a and a. To derive the 
equations of motion for these operators it will be conve- 
nient to have an expression for them in terms of the field 



operator 4'^. By normalizing 1$^'^) to unity we obtain 
a = j d'r<^T{r )i>a{r) (13) 

and 

5*,(r) = / dV'(f |Q-|r-")'I'a(r-"), (W) 



where Q'^^ is the projector on the subspace perpendicular 
to \^''^): Q^^^l - |$^'=)($^''|. To simplify the notation 
we introduce the operand o to denote the action of an 
operator O onto a field operator such as ^a'- 



cPrO\r)^a{r) 



(15) 



(o is also used in Ref. ||T^). With this notation Eqs. ( [T^ ) 
and (nj) can be written in the simpler form 



(16) 



Our calculation assumes a condensate of weakly inter- 
acting atoms {a^p <?C 1, where p is the atomic density) at 
very low temperature T « 0. In this limit the fraction of 
non-condensed atoms is very small, SN/N ~ 10^^ — lO^'^ 
p6| . The field operator has matrix elements scaling as 
\/ Na from the condensed atoms and contributions of or- 
der \/5Na from the non-condensed atoms, where Na and 
SNa are the average number of particles and the number 
of uncondensed particles of type a. Due to the small ra- 
tio between non-condensed and condensed particles it is 
useful to perform an expansion of the field operator '^a in 
terms of \/ SNa/Na- Formally this expansion is achieved 
by taking the limit 



NaQaa = COUStaut, 



(17) 



while keeping a constant ratio between the scattering 
lengths and between the number of particles in the a and 
the b states. Below we show that SNa is of order unity in 
this limit, and our expansion in terms of 5Na/Na be- 
comes an expansion in Na (note that since we have a 
fixed ratio Na/N), this also corresponds to an expansion 
in 1/V^). 



III. EXPANSION IN POWERS OF l/\/Na 

In this section we present the explicit expansion in 
powers of 1/ \/ Na ■ The expansion procedure is imple- 
mented by considering the time evolution of the operator 



A-(f,i) = 



-.a^{t)5i>a{r,t). 



(18) 
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From the definition of the condensate wavefunction (|1]^ 
and the splitting follows that the expectation value 
of A^'^ vanish exactly at all times 



(19) 



As we shall show below, A™ is of order unity in the 
limit in Eq. ( p7| ) so that the number of non-condensed 
atoms approaches a constant in this limit. To perform 
our expansion we write the operator A™ and the wave- 
function as a series of terms 




Na 



(20) 



and we derive equations for each of the terms order by 
order. Since we require that the exact wavefunction is 
normalized to all orders in 1/ \//Vq , the lowest order con- 
tribution must also be normalized ($a|^a) — 1- 

The procedure in this section is to calculate dA'^^/dt 
and da/dt to order k in Na (k — —1,0) and take 
the mean value (A^'^) which must vanish to all orders in 

1 / V Na- With this procedure we calculate the contribu- 
tions to the field operator to order N^. The main differ- 
ence compared to the derivation of Castin and Dum in 
Rcf. [p^ is that in their derivation the number of atoms 
of type a is fixed. This means that Noa — Na, where 
^Oa = a^a is the number operator for the condensed par- 
ticles, is only of order ~ 1 because the difference 
only arises due to the non-condensed particles. In our 
situation the distribution on different nu mbe r states is 
a binomial distribution with a width ^ \/ Na and hence 
the operator Noa — Na is of order \/ Na- 



In these expressions the only time derivatives are of the 
wavefunction $™ and the field operator ^a- By substi- 
tuting the splitting of the field operators in Eq. (^ into 
the right hand side of the time derivative of in Eq- (§) 
and inserting this into the above expressions we may find 
equations of motion for A"^ and d involving different pow- 
ers of 6^ a- Below we shall only consider the two lowest 
order contributions to A'^ and d, and it is sufficient to 
keep the terms which are linear in (5^a- Equivalently, 
this corresponds to neglecting the terms of higher than 
second order in the 6"^ operators in the Hamiltonian if 
we had used the approach of Gardiner where the 
approximation are made to the Hamiltonian rather than 
in the equations of motion. 



B. Order Na: Gross-Pitaevskii equations 

The largest terms in the expression in Eq. ( p^ ) are the 
th e te rms on the first line. These two terms are of order 
■\/ Na, but with a suitable definition of the wavefunction 
\^a^) the two contributions cancel each other to leading 
order in \/ Na, and A"'' is of order N^. By using Eqs. (^) 
and ^ and keeping only the dominant contribution from 
the condensed particles we get 



i^A- = ^/WaQa O (Uop^a - ) |$a) + 0{N"J, 



(23) 



where we have introduced the Gross-Pitaevskii Hamilto- 
nian 



a + 9aaNa\'i>a\'^ + 9abNb\<i>b\ 



(24) 



A. Basic equations for A^^ and a 

As a starting point for our calculation we use the time 
derivatives of the operators A'^ and d. By taking the 
time derivative of Eqs. (|l^) and (|l^) and using the defi- 
nition of A;^'^ ([l|) we obtain the expressions 



. d . 
I— a 
dt 



dt 



d ~ 



and 



dt " 



+( ( 1 oi<i>r^ 



dV 



08"^ a 



(21) 



(22) 



In Eq. (|2^) we have replaced the number operators by 
their leading order contributions a^a « Na and d)a^ dd « 
7V^. The difference between the leading order and the 
operators is of higher order and will be taken into account 
below. From Eq. (|2^) follows that the \fWa contribution 
to A';'' vanishes if we chose 



-Z- + HGP,a)|$a)=C(t)l<i>a) 



(25) 



The term on the right side introduces an optional phase 
factor on the wavefunction (C must be real due to 
the normalization condition) . Since all measurable quan- 
tities are expressed in terms of the field operator , we 
can introduce any phase factor on if we include the 
opposite phase factor on a. Here shall choose the sim- 
plest possible phase evolution of a and we set C = so 
that ( p5| ) reduces to the usual time dependent Gross- 
Pitaevskii equation. With this choice of C the equation 
of motion for d cancels to leading order idd/dt = 0{N^). 
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C. Order N°: = 0. 



D. Order A^^: Bogoliubov approximation 



In this subsection we shall show that the lowest order 
correction to the wavefunction <I>i^^ vanishes exactly in 
the considered experimental situation. To show this, we 
consider the mean value of Eq. (E2|) to order iV° 



*^(Aa(r)> =VNa{r\QTo ( Ho,a + gaaNa\<fT\' 



dt 



)+Oil/VNa). 

(26) 



When calculating the mean value of the third term in 
Eq. dH) {d'^l/dt di/aj/V^a we need to consider mean 

values like gaa{da^a^S^a)/\/N^- These terms do not 
vanish and at first sight they might seem to contribute 
to the present order. But since (A^'^) vanishes, this term 

may be rewritten as gaaiiNoa - Na)d'^6-^a)/VNa, and 
since the parenthesis is only of order y/ Na this term does 
not contribute to the present order. 

Because of the exact relation (|l^) the mean value of 
Eq. ( p^ must vanish to all orders in a/ Na- The ^/Nq, 
contribution vanish due to the Gross-Pitaevskii equation 
(p5|). The next order is obtained by expanding and 
'1>1'' in Eq. M) and we get 



= QaC 



+ <?,,iV,(2|$„|2$W+ci>2ci>(i) 



(27) 



Due to the projection operator Qa this equation admits 
the inclusion of a an overall phase in ^'"a\ i.e., id^^a^ /dt 
may have a contribution C'^'^a. Again we prefer to have 
the simplest possible equation for a and chose (^-^^ = 
0. With this choice the equation for ^i^-* is linear and 
homogeneous so that ^i^"* vanishes if it vanishes at i = 
0+. In Ref. it is shown that = in thermal 
equilibrium as we have just before the resonant pulse 
(|l|) , and since the condensate wavefunction is exactly the 
same before and after the pulse, as described in Eq. (12), 
we have 



$W(t) = 



(28) 



With the result ( Pq ) we can now calculate the deriva- 
tive in Eq. ( ^ ) to order and we get 

=.9aa($a||$aP|*a)(iVoa ~ Na)a 

dt 

+ 5afc(*a||$fcn$a)(iVofc-iVb)a 

+ 5aaV^((*a||$aP O A, + o \<^a\''\<^a 
+ 5afcy^((*a|*b*aoAf,^ 



Ajo$*$b|$„) |a, 



(29) 



where we have replaced expressions like a''d''da — d'^aNa 
with their lowest order contributions {Noa — Na)Na- The 
right hand side is seen to be a Hermitian operator multi- 
plied by d. This does not contribute to the time deriva- 
tive of the number operator Noa, so that we have 

j^Noa = 0{N^J. (30) 

This relation is consistent with Aa being of order unity. 
The number of particles of type a is a conserved quantity 

and hence Noa + 3-^10 5^ a = Noa + AloAa + 0{l/ y/K) 
must be conserved. Since Aa is of order N^, we must have 
a relation like (^0|). The weak time dependence of the 
number of condensed particles indicate that essentially all 
particles stay in the condensate modes and this justifies 
the description of the system by single mode approxima- 
tions. Note, however that such single mode calculations 
must be performed with care. A calculation which com- 
pletely ignores the non-condensed particles only contains 
the first two lines in Eq. (|29|), but the remaining terms 
are of the same order as the first and should be kept in 
the calculation. As we shall see below these terms gives 
rise to the same interaction as in the single mode approx- 
imation but with a different coupling constant. 
For the As the equation of motion reads 





'Aa- 




'Aa- 


d 


Ab 


= Co 


Ab 


I— 
dt 


At 


At 




UtJ 




UtJ 



where the vectors d and f3 are given by 



l*al 



.9aaV^I*a) 

gab\/%\^b) 

-gaaVN'alK) 

'N'bin) J 



(31) 



(32) 



and 



/3= |<i>fc| 



for all t. 




(33) 
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and where the matrix C is 



HGP,a + QagaaNa\'^a?Qa Q a9ab^ N aNb^ a^Qb Q a9aaNa<i>lQ*a Q a9ab\/%M'^a'^bQl 

Qb9ab\/NaNb^l<^>bQa HgY> .b + Qb9bbNb\^b?Qb Qb9ab^NaNb^a^bQ*a Qb9bbNb^lQl 

-Ql9aaNa^aQa ~Ql9ab^ NaNb^l%Qb ~ H GF ,a - Ql9aaN a? Ql ~Ql9abVWbK'^bQl 

-Ql9abVNaNbKnQa -QlgbbN^^'' Qb ~Ql9abV^M^anQl -~HGP,b - QlgbbNbl'^bl^Ql 

(34) 



Here we have introduced the projector Q* which projects 
onto the space orthogonal to |$*): Q* = 1 - |<i>*)($^| 
(note that Ql ^ Q\ = Q^). The first part of Eq. (||) 
containing the operator L corresponds to the equation 
derived for a single component in Ref . . With a fixed 
number of particles in the a and the 6 states this is the 
only contribution to this order. In a situation where we 
cannot ignore the fluctuations in the particle numbers 
there are additional contributions due to the last two 
terms in Eq. ( |3l| ) . Physically, we can understand the ori- 
gin of these terms by considering a subspace with a given 
number of particles Ma and Mb (of type a and h respec- 
tively). In this subspace the Gross-Pitaevskii equations 
are not the best approximations to the evolution of the 
condensate mode because they contain the average values 
Na and iVj, rather than Ma and Mb- From the structure 
of the last two terms in Eq. (|3^) it is seen that these 
terms corrects for the "wrong" value of the particle num- 
bers. The approximate method used in Refs. ^ is 
designed capture this dependence on the particle num- 
bers and in section ^ wc show that the results obtained 
with the present method are very similar to the result 
obtained in Refs. g, 



IV. SOLVING THE EQUATIONS 

With the equations (|||) and (|l]) we have derived the 
equations of motion to the desired accuracy, and in this 
section we solve these equations. The solution is obtained 
by expanding the A operators on a suitable set of mode- 
functions at i = 0+. The time evolution of these modes 
is chosen such that the first term in Eq. (^ij) is auto- 
matically taken into account and we then treat the last 
two terms in Eq. (31). Finally the solution for the A 
operators is inserted into Eq. (|29|) which is then solved. 



A. Expansion at t = 0"*" 

Instead of having the operators A^ and A;, which have 
both a spatial part and an operator character, it is con- 
venient to expand A^ and Af, on a set of vectors which 
take the spatial dependence into account. In such an ex- 



pansion, the expansion coefficients become operators and 
we have 



Aa 

Af, 
Ai 
At 



E 

fe=i 



Cfc 



" \Uak) ' 




r Kk) 1 


\Ubk) 
\Vak) 


+ 4 


Kk) 
Kk) 


. \vbk) . 




- Kbk) - 



(35) 



The UaS and the VaS can be any functions perpendicular 
to $a and $* respectively, but the whole set must be 
chosen such that the entire space perpendicular to the 
condensates modes can be spanned by the vectors. A 
convenient choice of the functions is discussed in section 

The u and v functions play the same role as the mod- 
efunction in ordinary Bogoliubov theory and it is an ad- 
vantage to choose the standard normalization 

{UakKk') + {ubkKk') ~ {VakKk') - {vbkKk') = 

{v*akKk') + (KkKk') - iKkKk') - {u*tkKk') = 0. 

(36) 

With this normalization we may express the expansion 
operators Ck in terms of the A operators 

Ck = {Uak\ 0^a + {Ubk\ O Ab - {Vak\ ° ^1 - {vbk\ ° ^l, 

(37) 

and from the commutation relation [A^ , A^] = Qa we 
find that the c^s obey bosonic commutation relations 

[cfc,4'] = ^kk' and [ck,Ck'] = 0. 



B. Time evolution 



To evolve the expansion (^5|) we need to choose the 
time evolution of the modefunctions. A convenient choice 



IS 



dt 



' Kk) ' 




' Kk) ' 


Kk) 


= Co 


Kk) 


Kk) 




\Vak) 


. Kk) . 




. Kk) . 



(38) 



With this choice of time evolution, the normalization ( |3q ) 
is conserved and so is the commutation relation of the Ck 
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operators. Indeed, the time evolution of Cfc only comes where the functions Fak and Fbk are defined by 
from the fluctuations in the number of particles and is 
given by 



ij^Ck = (Noa - No) fak + (iVob - Ni,)hk, (39) 



where 



Fak{t) = I dt'fakif) 



fak^{{Uak\,{ubk\,-{Vak\,-{vbk\) ■ a (42) 

Ibk={{Uakl{ubk\.-{Vakl-{vbk\)-l3. F^kit) = dt' hk{t') . 

J 

This equation is easily solved and by using that iVoa and 
N^b do not depend on time to this order of approximation 
we find 

Ck{t)^Ck{t^Q+)-i{Noa-Na)Fak{t) 

(41) Finally, by inserting Eqs. ^ and (gl|) into Eq. (11) we 
-i{Nob~ Nb)Fbk{t), obtain 



i-a ^{Noa " Na) (^gaa(<i>a||<faP|$a) + 2 Iml/^feF^fe) j a 



(43) 



{Nob - Nb) ( gabi^alM'fa) + 2 ^ Im(/:,F,fe ) j « + I E f:kit)ck{t - 0+) + fak{t)ci{t = 0+) ) a. 



k=l 



\k=l 



The expression in Eq. ([43|) has the solution 



a(t) = exp -z(77a(i) + Qa{t)) 



d{t = 0+), 



(44) 



where the Hermitian operator Qa{i) is just the time in- 
tegral from to <: of the three terms multiplying a on 
right hand side of Eq. (|43|). The phase 77^ arises because 
the term involving the Ck and c\ operators does not com- 
mute with itself at different times. By differentiating the 
expression in ( |4^ with respect to time and using the 
Baker-Hausdorff relation a few times, we find that (^) 
is indeed a solution of (Esh if the phase is given by 



Va{i)=Y. f dt'lm{Fak{t')f:k{t')). 

k=l -^0 



(45) 



A derivation of this phase and a geometrical interpreta- 
tion of it can be found in |2^. With the expression (^) 
we have finished our derivation of the evolution of the 
field operator. In the next section we apply the developed 
theory to describe phase collapse and entanglement. 



APPLICATION OF THE THEORY 



In this section we apply the theory developed in sec- 



phase between the two condensates and the entanglement 
of the atoms. For simplicity we shall not be as general as 
in the previous sections and we only consider a symmetric 
interaction which we describe below. 



A. Symmetric interaction 

We consider the symmetric situation where gaa — 
9bb 7^ ffob, and where the first pulse (|l|) is a tt/2 pulse 
such that Na = Nb- We assume a system with a fixed 
total number of particles A'' = 2Na so that we have 

Noa -Na= -{Nob - Nb) = (Noa " Nob)/2 + OiN^). 



tions II] and IV to describe the collapse of the relative 



(46) 



We also assume that the trapping potentials for the two 
different internal states are identical, spherically symmet- 
ric, and harmonic Va — Vb = . With this sym- 
metric choice of interactions the condensate wavcfunc- 
tions are exactly the same for the two components of the 
condensate $0 = = $. Furthermore the symmetry 
can also be exploited in the Bogoliubov modes. We di- 
vide the sum over k in Eq. ( p5[ ) into a sum over terms 
which are even under the exchange of a and b {+ modes) 
and terms which are odd under exchange of a and b (— 
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modes), i.e., we have 













" \Uak) ' 




' \<) ' 




" \nk) ' 


\ubk) 




H) 


or 


-\n-k) 


\Vak) 




H) 






. \vbk) . 




. H) . 







(47) 



In the remainder of this article superscript + and — on 
operators and functions will refer to these even and odd 
modes. The symmetry of the modes is reflected in the 
functions defined in section IV which obey 



fbk ~ ^ fak 



F, 



bk 



(48) 



From these relations we find that the two components 
have the same phase 



Va = Vb 



(49) 



and from the second relation and Eq. (46) follows that 
the operators describing the + modes are independent of 
time 



dt " 



B. Initial conditions 



(50) 



To evaluate the time evolution of physical quantities we 
need to relate the different operators after the resonant 
pulse (^ to the similar operators before the pulse. The 
relation between the operators before and after the pulse 
is completely describe by Eq. (|^) and below we extract 
some result from the general relations. 

Before the resonant pulse (|^) we assume that a con- 
densate is formed in the a state and that there are no 
particles in the b state. This means that d{t = 0~) is of 
order a/Iv whereas b(t = 0~) is of order unity so that we 
have from Eqs. (§) and (|l|) !§ 



a{t - 0+) - ^(a(O-) - &(0-)) « a{0-)/V2 
b{t = 0+) = ^(a(O-) + &(0-)) « a(0-)/V2. 



(51) 



These approximations should be handled with care since 
they violate both the commutation relation and unitarity. 
The sign « in Eq. (|5l|) means that the operators have 
approximately the same matrix elements. 

The operators can be related to the operators at 
i = 0~ by using Eq. (p7|). Inserting the approximate 
relations in (pTh we obtain 



(52) 



(i = 0+) = V2((4(0+)|oA„(0-) 

-(t;+(0+)|oAt(0-) 



where we have introduced a A operator before the 7r/2 
pulse Aq(0~) = a'l'(0")/ViV'5*a(0") which is also used 
in Ref. As the other A operators, this operator 

can be expanded on a set of Bogoliubov modes Aa(O^) = 
Lfc Cafe (0" ) (0~ )) +cjjj. (0~ ) l-Dafe (0" )) , whcrc the mode 
functions obey orthogonality relations similar to Eq. (|3^). 
In |]l^ it is shown that by a suitable choice of the 
modes |uafc(0~)) and |wafc(0~)) the pseudo-particle op- 
erators Cafc(0~) are the Bogoliubov operators describ- 
ing excitations of the condensates, and the ground state 
is the vacuum of these operators. The simplest rela- 
tion between the operators before and after the pulse 
is obtained by choosing |m^(0"'")) = |wafe(0~))/\/2 and 
I'^^fc (0"^)) = bafc(0"))/\/2. With these modes the opera- 
tors describing quasi-particle excitations in the + modes 
are exactly the same as the operators describing excita- 
tions before the pulse 



cj(i = 0+) =Cafe(< = 0-). 



(53) 



To find the initial condition for the operators we 
again use relation ( ^ ) and the approximation in Eq. ( |5l| ) 
and we find [Eol 



2fc(i = 0+) = - 



1 



a\^'){ul{^^)\obi>b{^-) 

-a(0-)(^;-(0+)|o5*t(o-)y 

(54) 

The simplest initial condition is obtained by choosing 
|wjr(0+)) = 0. With this choice c^(0+) is proportional to 
(5'I'6(0~). Since the initial state is the vacuum of (5^h(0^) 
the initial state for is also the vacuum state indepen- 
dent of the modefunctions |m^(0+)) and the state of the 
condensate. In our numerical work described below, we 
choose the ujTs to be the excited states of the Gross- 
Pitaevskii Hamiltonian before the pulse. 



C. Phase collapse 

We now investigate the collapse of the relative phase 
between the two components of the condensate. As a 
quantitative measure of the collapse we use the reduc- 
tion of the mean spin in the a;-direction which is defined 
by (Ja;) — Re((J+)), where J+ is given in terms of the 
fields operators by Jj^ = ^\o ^b- By inserting the split- 
ting in Eq. (^) , using the approximation in Eq. (|^) , and 
introducing the As we obtain 

(J+) = (fltfe) + (At o kb) + 0{\|^fW,). (55) 

The first term in this expression represents the time evo- 
lution of the condensed atoms and the second term is the 
contribution from the non-condensed atoms. At short 
times the first term is of order N„ whereas the second 
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term is of order unity. However, as we shall see below, the 
time scale of the collapse scales as ^/iVa, and since dka/dt 
is of order unity, the second contribution may become 
comparable to the first at the collapse time. We shall 
only work in the dynamically stable region gab < gaa- In 
the opposite case gab > gaa it is energetically favourable 
for the two components to separate so that any asym- 
metry in the two wavefunctions can grow with time and 
the system is dynamically unstable ||lj]. By using the 
same techniques as was used to show a similar result in 
[^5| , one can show that equation of motion for the us and 
vs (|3^) corresponds to the time evolution of a perturba- 
tion of the condensate wavefunction perpendicular to the 
wavefunction itself. Hence the stability of the system is 
determined by the stability of the solution of the Gross- 
Pitaevskii equation, and if the system is dynamically un- 
stable the number of non-condensed particles is expected 
to grow exponentially with time, so that the approxima- 
tions performed here are invalid. On the other hand, in 
the region gab < gaa it is energetically favourable for the 
two components to be overlapping so that the system 
is dynamically stable and the number of non-condensed 
particles performs small oscillations in time. This be- 
haviour is confirmed by direct numerical integrations of 
the equations. In the stable region the contribution from 
non-condensed particles remain much smaller than the 
contribution of the condensed particles, and in the fol- 
lowing we ignore the second term in E g. ( p5|). 

When we insert the solution in Eq. (^ ) into Eq. (|55| ) 
we see that the phase collapse arises from the difference 
of the angles Qa — ^t- By using the relations (ji^ ) and 
(Eq), the difference in the angles may be expressed as 



(56) 



where the function A is determined by 



- {gaa 



Jab){mm+^iMj2fak*Fak] (57) 



Kk=l 



and X{t = 0+) = 0. The first term in the time evolu- 
tion in Eq. ( pq ) is the same as the one obtained from 
the time evolution with a Hamiltonian like (^ . The cou- 
pling strength (|^) was derived by only considering the 
contribution from the condensate modes, and this cou- 
pling constant is the same as the first term in the correct 
coupling strength in Eq. (^7|). The second term in Eq. 
( p7| ) represents a modification of the coupling constant 
due to the Bogoliubov modes. The atoms in the con- 
densate mode does give the largest contribution to the 
Hamiltonian (H), but we do not obtain the correct cou- 
pling constant by only considering the terms from the 
condensate modes, because it is the fiuctuations of these 
terms, and not the terms them self, which give rise to 
phase collapse and squeezing. The fluctuations are of 
lower order and are comparable to the contributions from 




80 100 



FIG. 1: Modification of the coupling constant due to the Bo- 
goliubov modes. The figure shows the ratio 7 between the 
average of the correct couphng constant dX/dt ( |57[ ) and the 
average of the single mode value (^) for as^ab/o-s,aa = (+), 
as,ab/aa,aa = 0.5 (x), and as,ab/as,aa = 0.93 (A). The last 
value corresponds to the value for \F = 1, Mf — ±1) sodium 
atoms In the limit of very weak interactions {jj. ~ 3/20;), 
7 approaches the estimate in Eq. (^, and in the Thomas- 
Fermi limit {/J, 3> o;) it approaches 7/10 as predicted in Ref. 

0- 



the Bogoliubov modes, and the Bogoliubov modes have 
to be included. The Hamiltonian (^) is, however, still 
valid provided that we use the correct coupling constant 
in Eq. (13). 

The presence of the second term in Eq. (IstI ) only 
changes the coupling constant slightly. Rcf. ]14| consid- 
ers the same situation (by a different method) and in the 
Thomas- Fermi limit (^/w 3> 1, where fi is the chemical 
potential) it is shown that the effective coupling constant 
is approximately 7/10 of the coupling constant in Eq. 
(^. The functions and -FJ^ appearing in the second 
term in Eq. ( p7[ ) comes from the change in the condensate 
wavefunctions due to the fiuctuations in the number of 
particles in states a and b. Outside the the Thomas-Fermi 
regime the condensate wavefunctions are less affected by 
the interactions and are closer to the ground state of the 
harmonic potential which is independent of the number 
of particles of type a and b. We therefore expect that 
the contribution of the second term becomes smaller if 
we go away from the Thomas-Fermi limit. This is indeed 
confirmed by a numerical integration of the equations. In 
Fig. |l| we show the ratio 7 between the slopes of linear 
approximations to the time integrals of Eqs. (Et]) and (|^). 
In the non-interacting limit {n « 3w/2), Eq. (0) gives the 
correct coupling constant and in the Thomas-Fermi limit 
(/i ^ lli) the correct coupling constant (|5^) approaches 
7/10 of Eq. (^. Note, that this is only true for a time av- 
erage of the coupling strength. The expression in Eq. ( |57| ) 
have larger oscillation in time. 

If we ignore the non-condensed particles, the initial 
state is a Fock state with A'' particles of type a in the 
condensate mode and the vacuum of b and for all k. 



By using the exact relations in Eq. (51) we can calculate 
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the shape of the phase collapse 



N 



cos^ ^(A)exp 



k=l J 



(58) 



An example of the phase collapse is shown in Fig. ||, 
where we show the phase collapse obtained by numeri- 
cally integrating the equations. The functions F~^^ per- 
form small oscillations in time and the last exponential 
in Eq. ( p8| ) is of minor importance. The phase collapse 
arises from the cosine in Eq. ([s^). The value of (Jx) is 
reduced by a factor 1/e when A ~ a/2/ (N — 1) and the 
collapse time is roughly given by 



1 



7V-l7X' 



(59) 



where x is the average of the coupling constant in Eq. (|^). 
Since x is proportional to 1/A^, scales as VlV- 
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FIG. 2: Collapse of the relative phase between the two com- 
ponents of the condensate. The mean value of the spin 
(Ja:) decreases due to the fluctuation in the number of par- 
ticles in each of the internal states. The parameters are 
6 • 10"- 



N = 10^ 



as,bb 



2as 



and 



possible to optimize the measurement output by choosing 
the angle 6 such that the variance is minimized, and by 
minimizing the variance (AJg)'^ with respect to we find 
that the minimum variance is given by 



(AJ6 



)'• 

/mm 



and that the optimal angle is determined by 

i^J zJy JyJz^ 



tan(26l) 



(60) 



(61) 



To calculate the variance (AJ6i)miii it is an advan- 
tage to replace Jy by the raising and lowering operators 
Jy = (J+ - J-)/2i, so that, e.g., (J^) ^ (J+J-)/2 - 
Re(( J+ J+))/2, where we have used the symmetry to ob- 
tain (J^J_) = (J_ J+). By expressing the raising and 
lowering operators in terms of the field operators, u sing 
the splitting (^), and introducing the operator ( |T^ 
we obtain 

(J+ J_) = {a^a) + {a^ab^b) + Na{{Al o A,) + (A^ o A,)), 

(62) 

where we have neglected terms of order \/ Na- To cal- 
culate the first two terms in this expression we use the 
fact that we have conservation of the number of particles, 
and that Jz commutes with the Hamiltonian because we 
have no spin changing collisions. From the mean number 
of particles in the a state Na — {^a ° '^a) we find 



from ( J2) = N/4: = ((*t o - Na)^) we find 



(63) 



((ata)2) = 



N{N + 1) 



-N{AloAa)+0{VNa), (64) 



D. Spin squeezing 

We now turn to the calculation of spin squeezing. The 
calculation of spin squeezing is more complicated than 
the calculation of the phase collapse because squeezing is 
a lower order phenomena. The operators describing the 
noise are the square of the angular momentum operators 
and involve the expectation values of the product of four 
field operators. Such a product is in itself of order iV^ 
but due to interferences between different terms of or- 
der N'^ the noise is only of order N. Hence to calculate 
the squeezing we need to be very precise when we cal- 
culate the matrix elements. Here we calculate all matrix 
elements to order N and omit terms of order VlV and 
smaller. 

We are interested in estimating the optimal squeezing 
obtainable with the proposal B. Experimentally it is 



and with this expression and N^ — ((^pt o + 5*]^ o '^b)^) 
we get 



~N(AioK) + 0{y/Wa). (65) 



In the last equation we have made a replacement 
(di'lbbm'a) = N/2{Al o Aa) + 0{y/Na). This replace- 
ment is valid because 6 « a to lowest order according to 
Eq. (pll ). With these equations we can calculate (J+ J_) 
in (|62|) and by doing a similar calculation for (J4- J-|_) we 
find 



By similar arguments we can also calculate 



{JzJy 



JyJz) 



2Im( (a^atafe)) 



(66) 



(67) 
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The matrix elements appearing in Eqs. ( |6q ) and ( p% 
cannot be calculated to the desired accuracy by only us- 
ing Eq. (|^, because this equation does not go to high 
enough accuracy. We need to calculate the matrix ele- 
ments to order N but in Eq. (Q) we are omitting contri- 
butions to a of order 1/^/N. The next order correction 
to a gives a contribution of order N when inserted in 
Eqs. (66) and (|67|). However we can use the conservation 
of the number of particles in each of the internal states to 
find the remaining term in the matrix elements without 
going to higher order in the calculation of a. We split the 
exact operator a appearing in Eqs. (^6|) and (|67| ) into the 
part that we have calculated so far oq and an additional 
term Sa of order 1 / \/N 



ao 



6a. 



(68) 



At i = 0^ we have 6a = and by using that the number 
of particles in state a and aQOg are independent of time 
we find 



al{t)6ait) + Mt(i)ao(i) =At (0+) o A„.(0+) 

-At«)oA,(i). 



(69) 



Note, that we only need the lowest order contribution 
of this term, and since a{t) « a(0'*") « &(0"'") « b{t) 
according to Eqs. ( ^ ) and (^, we can actually replace 
So by bo which is necessary to complete the calculation 
below. 

With the above relation we are finally able to calcu- 
late the relevant matrix elements. Using the splitting in 
Eq. (|6^) , the initial condition in Eq. ( pl| ) , the time evolu- 
tion in Eq. (^4|), and replacing a{t = 0~)^a(t = 0~) with 
N - Aait - 0-yK{t ^0-)=N- 2A,(i - O+YKit = 
0+) we obtain 



2 _ NiN+l) 



N{N~1) 



cos^(2A)exp -8 5]|F-J 



(70) 



fe=i 



and 



{J^Jy + JyJ,) - cos^(A) sin(A) 




(71) 



By inserting these two expressions into Eq. ( pO| ) we can 
find the minimum variance and with the expression for 
(Jx) in Eq. ( p8[ ) we may find the squeezing parameter ^g. 

Note that we have not specified anything about the 
initial state of the -I- modes. Within the approximations 
used here, the results ( [70|) and (^ only involves ex- 
pectation values of operators for the — modes and the 
squeezing is independent of the state of the -I- modes. 
If the condensate is at a non-zero temperature before 




FIG. 3: Time Evolution of the squeezing calculated by a 
numerical integration (full line) and from the approximate 
Hamiltonian (g) (dashed curve). The parameters are the same 
as in Fig. Q. 



the pulse, it can be described by quasi-particle excita- 
tions, and these excitation are transfered into excitation 
in the -I- modes after the pulse. The — modes are al- 
ways in the vacuum state independent of temperature, 
and the squeezing calculated here is therefore also inde- 
pendent of temperature. Another nice property of Eqs. 
( [70| ) and ( |7l| ) is that they only involve the functions 
and not the number of non-condensed particles. Because 
the functions F^f^ vanish for states with non-zero orbital 
angular momentum, according to Eq. (^0|), it is sufficient 
to determine the evolution of the Bogoliubov modes with 
vanishing orbital angular momentum and this simplifies 
a numerical treatment of the system significantly. 

In Fig. 1^ we show the evolution of spin squeezing with 
the same parameters as in Fig. ^. In the figure we also 
show the prediction of the Hamiltonian (^) with a cou- 
pling constant x equal to the slope of a linear approx- 
imation to A. The two curves are in reasonably good 
agreement with each other. The deviation between the 
two curves is caused by the exponentials of the F~j^ func- 
tions in Eqs. ([t^ ) and ([n]), i.e., by the dependence of the 
wavefunction on the number of particles. After the res- 
onant pulse the wavefunction is no longer in equilibrium 
because the repulsion of the atoms is suddenly reduced, 
and the size of the atomic cloud will oscillate in time. 
The oscillations of the wavefunction depend on the num- 
ber of atoms in each of the internal state, and because 
the state of the system is distributed on states with dif- 
ferent number of atoms in each internal state, there is an 
uncertainty in the wavefunction which introduces noise 
and reduces the squeezing. After the completion of a full 
oscillation {ujt w 4, 9, 13, and 18), the wavefunction has 
approximately its initial form independent of the number 
of atoms in each internal state {F^f, « 0), and the results 
of the Hamiltonian (^) are in very good agreement with 
the numerical results. The result in Fig. |^ is very similar 
to the result obtained with the same parameters but a 
different approximate method in Ref. and our result 
thus supports the conclusions reached in that paper. 
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FIG. 4: Squeezing between the two \F — 1, Mf = ±1) hy- 
perfine states in sodium. The parameters are the same as in 
Fig. ^ except as,ab/ as,aa ~ 0.93. The full curve is the result 
of a numerical simulation and the dashed curve is from the 
Hamiltonian (^). 



The difference in scattering length used in Fig. || is 
an exaggeration of the reahstic parameters. If the ex- 
periment is performed with the two hyperfine states 
\F = l,Mp — ±1) in sodium representing the internal 
states a and b as proposed in Q, the ratio between the 
scattering lengths is as^b/as^aa ~ 0.93 |30 . In Fig. ^ 



we show the squeezing produced with this value of the 
ratio. Because the scattering lengths are very close, the 
oscillations of the wavefunction are much smaller, and 
the numerical curve and the result from Eq. (|^) are in 
better agreement. 

To create an entangled state the coupling constant x 
in the Hamiltonian ^ should be non-zero. So far we 
have obtained a non-zero coupling by having a difference 
in the the scattering lengths a^^aa + ^5,66 2as,a6- From 
the approximate coupling constant in Eq. (|^) we see that 
X is also non-zero if the a and the b atoms occupy two 
different regions in space, even if all scattering lengths 
are identical. In Rcf. |^ it was proposed to produce spin 
squeezing by replacing the two internal states a and b 
with two different momentum states which are created 
by a Bragg pulse shortly after the trap has been turned 
off, i.e., exactly the same setup as used in 31 . This 



proposal has the advantage that the squeezing is now be- 
tween atoms in the same internal state, and the squeez- 
ing is therefore better shielded from phase decoherence 
caused for instance by fluctuating magnetic fields. Also, 
this proposal could be used for rubidium atoms where the 
scheme which has been studied so far is not applicable 
because the scattering lengths are almost identical. 

To investigate this situation we assume that the time it 
takes the two momentum states to separate is very short 
so that we can neglect the interaction during this separa- 
tion process. In this case there is no interaction between 
the a and the b components and we can describe this sit- 
uation by settiiig a^^ab — and Va — Vt = after the 
pulse. In Fig. we show the result of such a simula- 
tion with parameters similar to the parameters in Ref. 



| pl[ . The calculation shows that strong squeezing can be 
produced with this proposal, but the agreement with the 
Hamiltonian (||) is not as good in this case because of the 
exponentials in Eqs. ( [70| ) and (|7l|). 



VI. CONCLUSION 

In this paper we have analysed the scheme proposed 
in Ref. |^| by describing the system with a number- 
conserving Bogoliubov theory. The developed theory is a 
consistent expansion in the ratio between non-condensed 
and condensed particles, and the validity of the calcula- 
tions can be investigated in a given experimental config- 
uration. 

The results obtained in this paper show that strong 
squeezing can indeed be produced by this method, but 
the theory is not able to determine the precise limit of 
the obtainable squeezing. The results of this paper agree 
with the result of the simplified Hamiltonian (^) which 
predict a reduction of the squeezing parameter by a fac- 
tor of approximately N~'^^^ , but we are ignoring terms of 
order a/ SN/N. From this we conclude that the obtain- 
able squeezing is at least of order ^J5N/N. The obtained 
results are independent of the temperature of the conden- 
sate, and the only effect of a non-vanishing temperature 
is that the approximations break down a little earlier. 
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FIG. 5: Squeezing created by a Bragg pulse as proposed in 
Ref. At time t = Q the trap is turned off and two com- 
ponents are separated by a Bragg pulse. The full line is 
the result of a numerical integration and the dashed curve 
is the prediction from the Hamiltonian (^. The parame- 
ters are as,aa/ \/ h/ M /uj = 6 • 10"", as,bb = as,aa, as.ab = 0, 
and — 1.7 • lO''. These parameters are chosen such that 
they are similar to the parameters in Ref. The ratio 

as,aa/ \/ h/M/cj = 6 ■ 10~* corresponds to sodium atoms in a 



trap with trapping frequency u = 2tv 
is turned off. 



19Hz before the trap 
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